Low lying twisting and acoustic modes of a rotating Bose-Einstein condensate 
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We present a calculation of the low lying spectrum of a rotating Bose-Einstein condensate. We 
show that in a cylindrical geometry, there exist two linear branches, one associated with usual 
acoustic excitations, the other corresponding to a twisting mode of the vortex lattice. Using a 
hydrodynamical approach we derive the elasticity coefficient of the vortex lattice and calculate the 
spectrum of condensate in a three dimensional harmonic trap with cylindrical symmetry. 
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By contrast with classical hydrodynamics, a quantum 
fluid cannot rotate like a solid body. Instead, it car- 
ries quantized vortices self organizing along a triangu- 
lar Abrikosov lattice when the rotation is fast enough. 
Vortices constitute a universal characteristic of quantum 
fluids and were directly observed in systems as different 
as liquid helium 0, superconductors Q, gaseous Bose- 
Einstein condensates (BEC) HHHH 

and very recently 

fermionic superfiuids 7] . The study of the excitations of 
a vortex lattice was initiated by the work of Tkachenko 
Q that was restricted to the study of the propagation 
of waves transversely to the rotation axis in an incom- 
pressible superfluid. The full understanding of these ex- 
citations in the case of dilute gases remains a challenge, 
due to the non trivial interplay between the elasticity of 
the vortex lattice and the phonon modes associated with 
the compressibility of these systems. Partial results were 
obtained in the case of a sin gle v ortex line 0, , two- 
dimensional systems 0, Hi! TlR| . fast rotating systems 
IIH E| or using the rotational hydrodynamics for- 
malism |l7l Hil Eof , but a complete theory remains to be 
found. 




In this letter, we present a study of the low lying modes 
propagating along the axis of an elongated vortex lattice. 
Up to now, these modes have only been considered for 
homogeneous and unbounded fluids [lf3 |. In the more 
realistic case of a trapped gaseous BEC, we show that 
these modes can be understood as Goldstone modes aris- 
ing from U(l) and 0(2) broken symmetries, associated 
respectively with the choice of the phase of the macro- 
scopic wave function of the condensate and of the direc- 
tion of the vortex lattice. These two broken symmetries 
thus give rise to two low energy excitation branches as- 
sociated respectively to a modulation of the phase of the 
wave function and of the direction of the vortex lattice 
in the z direction. Since the gradient of the phase is 
proportional to the velocity of the BEC, the hrst branch 
is associated with an acoustic wave propagating along 
the cloud ("phonon" branch). The second mode corre- 
sponds to the twisting of the vortices and is related to 
the elasticity of the lattice ("twiston" branch, Fig. fj). 
This scheme will be worked out first using a perturba- 
tive resolution of the Bogoliubov-de Gennes equations 
in a cylindrical trap, a method already used in the sta- 
bility analysis of solitons [l^. We will show that this 
solution can be interpreted in a macroscopic framework, 
leading to elasto-hydrodynamical equations for the mo- 
tion of the BEC in the presence of the vortex lattice. In 
particular, we propose the first derivation of the elastic 
response coefficient of the vortex lattice, starting from 
first principles. Finally, local density approximation will 
allow for the calculation of the lowest energy modes in 
the presence of an axial trapping. 

We consider a dilute Bose-Einstein condensate rotating 
at an angular velocity f^o along the z axis. At rest, the 
system is described in the mean-field approximation by 
a macroscopic wave-function i/jq solution of the Gross- 
Pitaevskii equation 



FIG. 1: Left: Structure of the twiston mode for a cylindri- 
cal BEC carrying three vortex lines. Right: arrangement of 
the vortex lattice in the planes defined by the shaded areas. 
The vortex lattice is rotated by an angle 8(z) depending on 
the axial position along the trap. This twisting induces the 
propagation of an elastic wave along the condensate axis. 



po + ff |^ | 2 )^o = 0. (1) 

where the single particle hamiltonian ho is given by 
ho = ~ ^-V 2 + V(x, y) — Q,qL z — fj, . V is the transverse 
trapping potential, L z the angular momentum in the z 
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direction, /xo the chemical potential and g the coupling 
constant characterizing the 2-body interactions. We as- 
sume for the moment there is no confining potential in 
the z direction. Nevertheless, we impose periodic bound- 
ary conditions in this direction, with period L, to allow 
for the normalisation of the wave-function. 

In the linear regime, we write function tp of the con- 
densate is given by ip — ipo + Sip. Sip is then solution of 
the Bogoliubov de-Gennes equations ihd t & = £<!>, with 
* = (Sip, Sip*) and 



C = 



ho + 2g\tp \ 



9^1 
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(2) 



We recall a few properties of C relevant for the dis- 
cussion num. First, although C is not hermitian, it 
is orthogonal for the quadratic form with signature (1,- 
1) defined by (4>i|<&2) = J d 3 r<fr\a3 i f>2, where CT3 is the 
diagonal Pauli matrix. Second, the eigenvectors $> a of 
C do not constitute a complete basis. Indeed, one can 
show that each zero energy mode is associated with an 
anomalous Jordan mode <b' a such that C&' a = & a . In the 
case under study here, we assume there are only two zero 
energy modes associated respectively with the U(l) and 
0(2) freedom of choice of the phase of the wave function 
and the orientation of the vortex lattice. The associated 
eigenvectors are respectively denoted by <& p (phase sym- 
metry) and <fr r (rotational symmetry) and are given by 

Mm 
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(3) 



Let us now proceed with the calculation of the long 
wavelength modes of C. Using the translational in- 
variance of the system we can write <f> a (x,y,z) — 
& a (x, y)e lkaZ . We therefore see that C = Cq + SC. 
where Cq acts on the transverse degrees of freedom only, 
SC = £fc a (T3, with ek a — S 2 i„/2m. 

Since we only care for long wavelength eigenmodes, 
we have ek a vanishingly small and we can treat SC as 
a perturbation of Cq. However, starting the perturba- 
tion expansion from the zero energy modes <fr r and <& p 
of Cq associated with the rotational and phase symme- 
tries gives rise to linear excited branches. Indeed, due 

to the non diagonalizability of £0 the first order term 

1/2 

of the low momentum expansion scales like e fc rather 
than £fc a [2^ . This singularity can be proven rigorously 
by diagonalizing the projection of SC on the vector space 
spanned by the zero energy and anomalous modes. Phys- 
ically, this result corresponds to the fact that one should 



recover a linear phonon dispersion law E a cx k a oc e]f 2 
associated with usual acoustic waves. 

Working out the perturbative expansion in power of 
y/tkl yields 

* a = (a r * r + a p ® p ) + E a (a r & r + a p * p ) + (4) 
where E a is solution of the eigenequation 



7— r,p 



( e fe a ^4/37 — E^Bp^dy — 0, 



(5) 



with (3 e {r,p}, Afa = (^^\a z ^^)/L, and Bp n = 
($p\&')/L. L is introduced here so that A and B are 
defined as ID linear quantities that we can use in a lo- 
cal density approach as demonstrated later in this paper. 
Using Eq. [21 yields the simple expressions 



(* p |a,* p ) = 2N (*„!*;) = d^N 

(* r |ff,* r ) = 2(1%) (* r |*;) = a f2o (£ 2 ) 

(* p |a,* r ) = 2(L Z ) (#„!#;) - 8n N 

(* r |a z * p ) = 2(£ z ) (* r |* p ) = a W) (L 2 ), 



(6) 



where N is the total atom number |2J]. The set of 
two equations JSJ yields two different linear excitation 
branches. The highest branch has a non-zero velocity 
for vanishing f^o and can therefore be identified with a 
phonon mode. The other branch has a vanishing velocity 
at small fl and is therefore the twiston mode. 

It is striking that the coefficients of Eq. (J5J can be 
expressed as simple combinations of macroscopic quanti- 
ties such as the atom number or the angular momentum. 
Similarly to superfluid hydrodynamics, this suggests that 
the formalism developed here can be expressed in term 
of an elasto-hydrodynamical theory. In order to clarify 
this link, we first note that the frequency of the long 
wavelength modes we are interested in is much smaller 
than the frequencies characterizing the evolution of the 
transverse degrees of freedom. This means in particular 
that the dynamics in the (x, y) plane is frozen and that 
we can therefore define local chemical potential fj,(z,t), 
angular velocity f2(z, t), phase of the wavefunction %(z, t) 
and angle of the vortex lattice 9(z, t). In other word, the 
wave function can be written as 



ip(x, y, z, t) = tpo(fi(z, t),fl(z, t),x(z, t),6(z, t),x, y) (7) 

where ipQ is the solution of the Gross Pitaevskii equation 
Q). Expanding the macroscopic quantities around their 
equilibrium values yields 



Sip = S/jd^ipo + Sfldn ipo + Sxd Xo ip Q + S0de o ip o , (8) 
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where /i = fiQ + Sfj,, fl = flo + S^j,, etc. Noting further that A similar analysis shows that the a — r component of 

5ip is the first component of 3> Q and that d x ipo = iipo and Eq. JSJ is equivalent to the angular momentum conser- 
d 9 ip = -iL z ip /h, we see that Eqn. © is equivalent to vation if the current T is taken to be 



* = i5x&p - i-r *r + 8(J,& p + 5n<f>' r . (9) 
Comparing to Eq. (@J this leads to the identification 

a r = -i80/h = Sn/E a (10) 
a p = iS\ = S^/E a , (11) 

which gives, using the identity dt — —iE a /Ti, 

d t se = sn (12) 

hd t Sx — ~S^l (13) 

The interpretation of these two relations is straight- 
forward. Indeed, ea. |12l is the analogue of the Kelvin- 
Helmholtz theorem |25| and implies that the vortex lat- 
tice rotates at the same speed as the local velocity flow. 
As for eq. El it is a version of the Bernoulli theorem, 
with hSx/m playing the role of the velocity potential. 

Starting from the identification of the a TtP coefficient, 
we now show that the eigensystem {^J can be interpreted 
as conservation laws for the system. Indeed, the conser- 
vation of particle number, momentum and angular mo- 
mentum in the z direction yield at first order in pertur- 
bation the very general set of equations 

d t 5n + d z {n o v) = (14) 
mn d t v + d z f = (15) 
d t S£ z + 8 Z T = (16) 

where uq is the linear particle density, v is the local ve- 
locity in the z direction, l z is the angular momentum per 
unit length. / and T are respectively the ID momentum 
and angular momentum currents, which can be identified 
with a force and a torque and will be calculated from the 
microscopic formalism presented above. 

Combining Eq. (|14|l and 1|15J) and using a generalized 
Gibbs-Duhem relation |26j for a rotating system df = 
nod/i + £ z dQ yields 

d 2 [Sfidp n + Sfldfi n ] = —d z [n d z 5fi + £° z d z Sil] 

(17) 

where we have used the adiabatic following of the trans- 
verse degrees of freedom to set 5n — S^id^nQ + 5SldQ o n . 

After making the replacement Tidt — —iE a and d z — 
ik a , we see that Eq. (|17f) is strictly equivalent to the 
a — p component of Eq. JSJl . 



r = i° z v - ndje, (is) 

where k = AL 2 Z /Lm, and AL 2 Z = (L 2 ) - {L z ) 2 /N is 
the fluctuation of the angular momentum. The term of 
r proportional to the local velocity v corresponds the 
convective part of the angular momentum current. The 
second term is associated with the elastic response of the 
vortex lattice to a torsion and k is therefore the elastic 
modulus of the lattice. 

The hydrodynamical formalism developed above can 
be used to extend our calculation to the case of a weak 
trapping in the z direction. In this case, local den- 
sity approximation can be applied. At equilibrium, the 
angular velocity is still uniform along the cloud. By 
contrast, the local chemical potential is now given by 
fj,o(z) = fi c — muj 2 z 2 /2, where \i c is the chemical poten- 
tial at the center of the trap. This permits to define local 
matrices A(z) and B(z) and yields the equation for the 
vector X(z, t) = (5^i(z, t),5Q(z, t)) 

2md 2 [B(z) ■ X(z, t)] = d z [A(z) ■ d z X(z, t)} (19) 

In general, A and B must be calculated using a numer- 
ical resolution of the Gross-Pitaevskii equation. How- 
ever, in the regime of fast rotation, quantities appearing 
in Eq. H19|) can be evaluated assuming a classical rota- 
tionnal flow v = J~2o x r. We have first checked that at 
zero angular velocity we could recover the spectrum of 
an elongated condensate studied in j^. We have also 
verified that the center of mass was indeed evolving at a 
frequency u> z , in accordance with Kohn's theorem. The 
other eigenfrequencies uj a are evaluated using a varia- 
tional method. We indeed note that the u) 2 a are the 
eigenvalues of the operator T = B^ 1 ■ d z (A ■ d z -)/2m 
which is hermitian for the scalar product defined by 
(Y\X) = J Y*(z) ■ B{z) ■ X(z)dz. We have calculated 
the ui a using polynomial trial wave functions of order 5. 
The variation of to a for the two lowest twiston modes - 
corresponding respectively to odd and even symmetries 
HI - is displayed in Fig. El 

In conclusion, we have demonstrated the existence of 
two low energy excitation branches of a rotating elon- 
gated Bosc-Einstein condensate. The validity of our cal- 
culation is limited by two conditions. First the applica- 
bility of the hydrodynamical formalism to evaluate the 
matrix elements of A and B requires a dense vortex lat- 
tice, hence a large rotation frequency. Using an imagi- 
nary time evolution of 2D Gross-Pitaevskii equation [29j , 
we have checked it was true for fio ^ 0.8w_l. Second 
the breakdown of the perturbation expansion will hap- 
pen when the tu a become comparable with the frequency 
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FIG. 2: Low lying twiston spectrum of a trapped rotating 
Bose-Einstein condensate. The two branches correspond to 
the frequency of the lowest odd (full line) and even (dashed 
line) twiston modes and were calculated using polynomials 
trial wave functions. 



lot of the first transverse excited mode, i.e. the low- 
est Tkachenko mode. Using the expression of lot found 
in 0], we find that for typical values n/hoj± = 20 
and lo±/lo z = 20, the perturbation expansion fails for 
rio/uj± < 0.96. Although the width of validity of our 
approximation might seem narrow at first sight, it must 
be noted that experimentally regular vortex lattices are 
only observed on a relatively narrow range of rotation fre- 
quencies, comparable to our validity domain (Qq/loj_ G 
[0.7,0.96] for [2fj). A last issue concerns the experimen- 
tal observation of twistons. It must be noted that, in 
practice, it is difficult to excite directly the vortex de- 
grees of freedom of a gaseous BEC and it is therefore 
more convenient to excite an acoustic mode by perturb- 
ing the trapping potential and to use the linear |3ll l32j 
or non linear |33j coupling to vortical modes. In princi- 
ple, the coupling existing between acoustic and torsional 
degrees of freedom should permit an excitation of the 
twiston mode using an external potential. However, one 
can show using perturbation theory that a simple mod- 
ulation of the trap center or frequency is only weakly 
coupled to the twiston mode. To obtain a significant cou- 
pling to torsional modes, one therefore needs to impart 
cubic or quartic perturbation to the gas. 
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